User guide of GranulatShiny
Authors : Aurel Hebert–Burggraeve, Laure Simplet, Laurent Dubroca,
Camille Vogel
Homepage
GranulatShiny is an application that facilitates the statistical
processing of data collected as part of the initial studies,
pre-construction baseline studies and environmental monitoring dedicated
to fishery resources and ichthyofauna relating to the appraisal of
applications for authorisation to extract marine aggregates. The
application automates some of the data formatting and the calculation of
standard biodiversity indicators, and provides decision keys for the
more advanced processing stages.
Based on the calculated indicators
and the user’s choices, the application produces figures and tables in
formats corresponding to the recommendations of the french reference
documents: the “Halieutic Protocol” and the methodological guide for the
development of Orientation Documents for Sustainable Management of
Marine Aggregates (DOGGM). The application provides an interactive
graphical interface based on the R language, relieving the user of the
need for mastery to focus on parameters of interest for diagnosing the
potential effects of marine aggregate extraction on fish resources.
GranulatShiny consists of 3 statistical approaches (exploratory,
descriptive, and inferential) that can be used primarily to quantify the
influence of marine aggregate extraction on fish communities. Each
approach covers an expected aspect of the “Halieutic Protocol”. The
exploratory approach presents and analyzes data at the scale of the
entire community. The descriptive approach presents and analyzes data at
the scale of a species. And resorting to inferential analysis is
necessary to evaluate the temporal and spatial variability of different
indicators of fish resources before and during exploitation.
The
current version of the application does not have a regulatory purpose
but serves as an aid in the production of monitoring reports on fish
communities. It does not replace the work already provided by consulting
firms but complements what is already done by allowing the assembly of a
statistical model to test the effects of certain parameters on fish
communities.
This guide has been written to enable a GranulatShiny user to
familiarize themselves with the application interface and understand the
methodology developed behind each result provided by the
application.
Before starting, note that there are different buttons
in the application.
The buttons with a boat icon allow you to switch from one tab to another.
Those bearing a small green dragon indicate a mandatory stopover.
The buttons with an arrow allow you to download results from the application. The formats used are (csv, png, txt, rds).
Finally, those with a circle containing an i are helps that can be displayed to help you better understand a graph or other objects proposed by the application.
When launching the application, the homepage opens automatically. On
this page, various reference documents and information are listed with
their associated URLs hyperlinked. A reminder of the context regarding
marine aggregate extraction is located on the right side of the page.
To move to the next tab, press the “start” button on the homepage
within the application.
Data formatting
Information to enter: loading datasets to analyze.
The first step in using the application is to import the data
collected. The data must comply with a standard defined by Ifremer, the
outline of which can be found here:
https://raw.githack.com/GranulatShiny/GranulatShiny/main/Description_Format_Generique_GranulatShiny.html
There are three possible scenarios:
- Case 1: You are new to the tool and have no formatted data
with which to test its functionality. In this case, you will use the
dataset supplied with the tool
If this is an exercise in
discovering the tool, a fictitious dataset is made available to the user
along with the tool. The data is directly integrated into the
application and can be loaded by selecting the “Non” answer under the
heading: “Do you have your own data? This dataset is intended for
educational purposes. It does not correspond to any real case and cannot
therefore appear in documents with administrative value (i.e. monitoring
reports, initial environmental status, reference status before works,
etc.).
For the purposes of familiarisation with the tool, the dataset made
available with the tool corresponds to a fictitious concession located
in the Bay of Biscay. For this fictitious dataset, we consider a
concession in operation from 2000 to 2030, for which monitoring of the
fisheries compartment has been set up every 5 years with 2 years of
initial status. The fictitious sampling plan provides for the sampling
of 10 stations within the concession and 10 stations outside the
concession. This choice does not correspond to a sampling plan that
would have been developed with knowledge of the environmental conditions
of the site (i.e. sedimentary facies, benthic habitats) and therefore
does not correspond to the recommendations of the “Halieutic Protocol”.
This fictitious sampling plan uses a beam trawl with a horizontal
opening of 4.4 m and a tow length of 1000 metres.
To avoid any
confusion, the 4 species present in this dataset are also fictitious.
Each species has a population dynamic associated with a specific
probability distribution law. Knowing this, it is possible to control
the results from inferential statistics and the effects of the
environment on the chosen species. The first species,
Cephalaspis.tenuicornis, was not affected by the extraction and its
spatio-temporal dynamics were stable throughout the monitoring period
(i.e. there was no effect of time, space or environmental conditions on
the observed abundances). Thus no potential effect of the variables on
the abundance of this species will be detected. The Dimichtys.terreli
species is impacted by extraction, but its spatio-temporal dynamics are
stable throughout the monitoring period. There is therefore a
significant effect of extraction, which is reflected in the difference
between the values obtained by sampling in or outside the exploitation
zone. Leedsischthys.problematicus is affected by extraction, but
differently depending on the season. The seasons do not influence the
population of this species in normal times, i.e. in the absence of
aggregate extraction, but the interaction between the effect of
extraction and the seasonal effect modifies the abundance of this
species. Finally, Latimeria.chalumnae is not affected by the effect of
extraction, nor by spatio-temporal effects such as season and/or
environmental conditions, but the abundance of the species is naturally
highly variable. These examples illustrate different responses to the
environment in order to better understand what is sought during
inferential analysis.
- Case 2: you have some initial experience of the tool and
are starting to analyse your own dataset in the recommended
format
If you are analysing real data in the appropriate
format, you need to select and load the following files into the
GranulatShiny graphical interface: “TuttiCatch.csv” and
“TuttiOperation.csv”, which contain most of the information relating to
the progress and results of the monitoring carried out. Only the csv
format is supported. The “TuttiCatch.csv” file corresponds to the catch
data from the sampling of fish populations and the “TuttiOperation.csv”
file corresponds to all the information derived from the implementation
of the protocol for each sampling station (i.e. date and name of the
survey, fishing gear, the characteristics of which will be specified in
the reports associated with the results, geographical coordinates of the
spinning and turning points and associated times, total duration of the
trawl haul, turning and spinning depth).
WARNING.
The expected data format must be respected, otherwise the processing
routines cannot run correctly and a warning message will appear on the
interface. In this case, it is recommended that you review your file
format with the expected file format. Furthermore, in the case of
sub-sampling or information provided at individual level, it is
important to report the information at tow scale so that each
combination of species and tow corresponds to a single line in the data
table. Otherwise a message will be returned indicating the existence of
duplicates preventing the catch data from being processed.
Once the files have been loaded, you will have access to a number of
new functions. A map centred on the concession will appear, displaying
the sampling stations. You will also be able to interact with the
‘Impact stations’ and ‘Reference stations’ fields. You will also be able
to import ‘ShapeFiles’ to display the contours of the marine aggregate
extraction concession.
Under the heading “Impact stations”, you can
check and modify the operating period (the period during which
extraction work takes place). You must also enter in the corresponding
space the stations that are impacted by the extraction. The colour of
the various samples will then change to red for the stations affected
(see figure below).
In accordance with paragraph 8.3.2 of the
“Halieutic Protocol”, the application is developed for the most common
case of trawl sampling. The horizontal opening length of the trawl must
be entered to calculate the sampled areas in order to work in density.
At this stage of development, the application does not take into account
other fishing gears.
Finally, in the exceptional case where a station
entered in the “TuttiOperation.csv” file needs to be removed after the
fact, this can be done under the “Reference stations” heading.
- Case 3: You have already used the tool to process your
data. You have a summary file of all the parameters used for a previous
analysis and you want to start again from this file.
If
you have already saved the settings in a file, you can import them after
the “TuttiCatch.csv” and “TuttiOperation.csv” files, so that the station
fields are filled in automatically.
Production of indicator tables
When you have completed the data loading stage, you can press the
button with the green dragon. This will launch the internal calculation
of the various indicators and covariates required to analyse the data.
If you do not press this button, nothing will happen and you will not be
able to continue with the analysis.
Nota Bene : If you have more than one concession to analyse, you can return to this tab, change the files by loading those corresponding to this other concession (“TuttiCatch.csv”, “TuttiOperation.csv”, “ShapeFiles”), then press the green dragon again to restart production of the indicator tables.
General table
In the “Tables” tab, there is a data table on the right and an interactive section on the left. The table displayed is formed from the data entered in the ‘Data formatting’ tab. The table formatting functions will calculate abundance, biomass and various diversity indicators for each station and for each survey. The “treatment” variable, which indicates the state of each station, can take two values: “no impact” or “impact”. It indicates whether the station is within the perimeter of the concession and therefore considered to be impacted by the aggregate extraction work (i.e. “impact” mode) or whether the station is outside the perimeter of the concession (i.e. “no impact” mode). In the case of an initial state, where there has been no extraction on the site of the concession studied, the stations located inside the concession are assigned the “no impact” state until the start date of exploitation. This allows them to be considered as reflecting the state of the environment before any impact from extraction, for the purposes of the statistical analysis carried out afterwards.
During the formatting processes, the season mode is calculated on the basis of the sampling start dates. The administrative framework is chosen by default to determine the seasons. However, it is possible to change this column in the full table. Particular attention is paid to the notion of season, as this is an integral part of the assessment of the temporal variability of fish communities (paragraph 8.1.4 of the Halieutic Protocol). According to the “Halieutic Protocol”, the effects of seasonal variability on fish assemblages (groups of species) are highly dependent on latitude. In northern waters (North Sea, English Channel, northern Bay of Biscay), it is common to observe only two types of fish assemblages per year, a winter assemblage for about eight months of the year and a summer assemblage for about four months. In the warmer waters of the south (south of the Bay of Biscay, Mediterranean), seasonal assemblages are potentially more numerous, with more marked spring and autumn assemblages. Nevertheless, it is the conduct of the initial survey that will make it possible to determine the seasonal variability locally and to decide on the seasonal periodicity of the monitoring surveys. By modifying the “season” column in the table, it is possible to adjust to local conditions and sampling difficulties.
You can change the general display of the data table using the arrow under the “which table to display” message, and you can download the table displayed using the “Download table” button. The “Download entered information” button is used to save the list of impact stations, the exploitation dates and the trawl opening width used in the “Import data” tab in a csv file. At the end of this tab, you can decide whether to carry out the “exploratory statistics” section, which looks at the community as a whole, or to go straight to the “descriptive statistics” section, which focuses on a specific variable.
Exploratory statistics
Representation of indicators
This section looks at the biodiversity and abundance indicators obtained at the community level, inside and outside the concession, for each data collection survey carried out. The indicators presented are those referred to in the “Halieutic Protocol”, article 8.4.1.
Definition of indicators
Biodiversity encompasses the variety of life at all levels of organisation, classified according to evolutionary (phylogenetic) and ecological (functional) criteria. At the level of biological populations, genetic variation between individual organisms and between lineages contributes to biodiversity as a signature of evolutionary and ecological history and a basis for future adaptive evolution. It is at the species level that the term biodiversity is most often applied by ecologists and conservation biologists. Species richness refers to the total number of species present in a given ecosystem. It is a simple measure that only takes into account the number of species without considering their relative abundance. For example, if a tropical forest contains 100 different tree species, its species richness would be 100.
A diversity index is a mathematical expression that combines species richness and evenness to measure diversity. The main objective of a diversity index is to obtain a quantitative estimate of biological variability that can be used to compare biological entities in space or time. This index takes into account two different aspects that contribute to the concept of diversity in a community: species richness and homogeneity.
The Shannon-Weaver diversity index is a widely used index for comparing diversity between different habitats. It assumes that individuals are randomly sampled from a large independent population and that all species are represented in the sample. This index measures both species richness and the equity (or uniformity) of species distribution in an ecosystem. It takes into account both the number of species present and their relative abundance. More specifically, the Shannon-Weiner index is calculated using the following formula:
\[ H′= -\sum_{i}^S (p_{i}*ln(p_{i})) \]
où :
S is the total number of species,
pi is the proportion of the i-th species among all the
species present,
ln represents the neperian
logarithm
The value of the Shannon-Weaver diversity index is generally between 1.5 and 3.5 and rarely exceeds 4.5. A higher Shannon-Weiner index indicates a greater diversity of species and a more uniform distribution between these species.
Unlike the Shannon-Weiner index, the Simpson index focuses primarily on the dominance of the most abundant species in an ecosystem. It is calculated using the following formula:
\[ D= \sum_{i}^S (p_{i}*(p_{i}-1)) \]
where the terms are the same as in the Shannon-Weiner index. A higher Simpson index indicates lower biodiversity, because it places greater emphasis on the probability that a species chosen at random is the same as the one chosen previously.
Representation in the application
Firstly, the table (below) shows the mean values for abundance, biomass, species richness, Shannon and Simpson indicators inside the concession, outside the concession and overall for each survey. They are calculated from the values obtained at each sampling station. For easier reading, the standard deviations are not displayed in the table in the application but are available in the cvs file that can be downloaded via the “Download table” button.
The graphs below show the mean values (dots) and the 5 and 95 percentiles (high and low bars) obtained for the same indicators as those in the table, depending on the survey selected and the sector sampled (paragraph 8.4.1 of the Halieutic Protocol). They provide a quick overview of the differences in values obtained between the concession area and the reference area for the most common biodiversity indicators.
The advantage of these approaches is that the fish community can
be compared on several scales. Initially, the comparison focuses on the
inside or outside of the concession. But if the surveys are looked at
one after the other, it may be possible to distinguish changes over
time. There is both a spatial and a temporal aspect.
Representation of the structure
This table represents the proportion of each species present for each sampling survey (paragraph 8.4.1 of the Halieutic Protocol). The table makes it possible to monitor changes in the proportions of species over time and provides a perspective on trends in assemblages. Variations in the proportions of different species from one year to the next can indicate significant ecological changes, such as fluctuations in biodiversity, changes in habitats or environmental pressures. This table can also be used to identify which species are dominant in a given ecosystem and which are in decline.
The figure below is made up of three graphs representing the abundance of species in a given survey in descending order (left), the relative contribution of each species to the total abundance (top right) and the species accumulation curve (bottom right). It provides information to meet the expectations of paragraph 8.4.1 of the ” Halieutic Protocol “.
The figure on the left shows the abundance of each species for a
survey in a histogram ordered by decreasing abundance. This makes it
possible to identify the dominant species within a survey, i.e. those
that are most abundant in the sample. These dominant species can play a
crucial role in structuring the ecosystem studied, influencing, for
example, competition for resources or predation on other species.
In addition, by observing changes from one survey to the next, this
histogram can help to visualise temporal trends. For example, an
increase or decrease in the abundance of a dominant species could
indicate changes in environmental conditions or in the interactions
between species. It can also be used to detect seasonal and annual
variations, such as breeding peaks or seasonal migrations, which can
influence the composition of the community.
Comparing data between
different surveys is also made easier by this histogram. By placing
species abundance distributions for different periods side by side, it
can help to identify similarities and differences between ecosystems at
different times of the year. This comparison can reveal general
ecological patterns or specific responses to environmental disturbances.
For ease of reading, species accounting for less than 1% of total
abundance are not shown on the histogram.
The figure in the top right corner is a curve of cumulative abundance as a function of the number of species in descending order of abundance. Cumulative abundance refers to the cumulative sum of species abundances in a dataset, starting with the most abundant species and successively adding the abundances of the following species in descending order. The cumulative abundance curve allows for the assessment of species diversity and distribution in an ecosystem or biological sample. The straighter the curve, the more diversified the community, whereas a curve that rises quickly then flattens indicates a community where a few species are very abundant while most species are rare. It’s important to note that the cumulative abundance curve is of interest when there are many different species. Here, the curve is constructed from the dataset of a fictitious concession with only 4 species. Therefore, it has limited relevance. In practice, you should not have this kind of result with your data.
The figure in the bottom right corner represents the number of
species as a function of the number of sampled sites. The algorithm, for
a given number of sites, will test all existing combinations in the
dataset and retrieve the number of species for each combination. Then it
calculates an average per number of sites, and it’s this average value
that is plotted on the graph.
The shape of the species accumulation
curve can provide information about the diversity of the studied
ecosystem. If the curve increases rapidly and tends towards a flat
asymptote, it suggests that most of the present species have been
sampled, providing a good estimate of ecosystem diversity. Conversely,
if the curve increases slowly and does not seem to reach a plateau, it
indicates that there are still species to be discovered, and sampling
should be continued to obtain a more accurate estimate of diversity.
The species accumulation curve is useful for determining the minimum
number of samples needed to obtain an adequate representation of species
diversity. This can be used during the initial assessment phase to
ensure that the proposed sampling plan captures the diversity of the
area. Here, constructed from the fictitious dataset, the curve holds
little interest.
Descriptive statistics
Data representation
In this section, we focus on a specific indicator (species abundance,
total biomass, diversity index, etc.) and compare it to the explanatory
variables in our dataset. We are looking for potential effects or
correlations upstream of inferential statistics. Therefore, we need to
select a variable that we aim to explain based on parameters related to
data acquisition and extraction. The statistical analysis will be
conducted on this variable.
Initially, the table summarizing the
explained variable provides information on the number of zeros and
missing values, the total length of the value series, and the fraction
of zeros and missing values compared to the total values. Additionally,
it also gives the mean, extremes, standard deviation, and quartiles of
the series. The most important values for the “inferential statistics”
part are those indicating the number of zeros and missing values
(“n_missing”, “complete_rate”). The proportion of zeros in the data
structures the modeling method to be employed. Indeed, a high proportion
of zeros can compromise the implementation of a generalized linear mixed
model.
Following this table, it’s possible to visualize boxplots. The boxplot provides another representation to interpret the relationship between the explained variable and the explanatory variables such as impact, year, campaign, station, and season. In graphical representations of statistical data, the boxplot is a quick way to illustrate the essential profile of a quantitative statistical series. The boxplot summarizes some position indicators of the studied characteristic (median, quartiles, minimum, maximum, or deciles). It is often used to quickly compare two series. In GranulatShiny, the series of the explained variable (here abundance) in the impacted zone is compared with that of the non-impacted zone. In the figure below, the same variable is represented on a decimal scale (on the left) and on a logarithmic scale (on the right). The logarithmic scale is offered in the application to transform the variable into a pseudo-normal distribution and provide more meaning to the boxplot representation by limiting the influence of extreme values. Each boxplot is constructed as follows: the horizontal line crossing the white square corresponds to the median, the upper and lower edges of the white square correspond to the 75th and 25th quartiles respectively, the ends of the “whiskers” correspond to the 95th and 5th percentiles; finally, the points correspond to extreme values.
In paragraph 8.4.1 of the “Halieutic Protocol”, it is indicated that the data should be described and analysed by size group, maturity or functional group. The GranulatShiny application does not allow data tables to be subdivided into subgroups. This must be done upstream by the user. To analyse a particular functional group, you need to sort your “TuttiCatch” file and save this new table so that you can integrate it into the application. This enables the GranulatShiny statistical method to be applied to the species to be treated separately.
Once you have explored the data, you can move on to the next tab by pressing the “Choose distribution probability” button or by clicking on “Diagnosis of analysis”.
Diagnostic analysis
This tab allows you to select and visualize the probability
distribution that best fits the explained variable. The frequency
histogram of the variable under study (selected by the user) is
represented by gray bars. It depicts the empirical distribution of the
observed data. It’s constructed by grouping the data of the variable
into intervals and counting the number of observations in each interval.
This provides a visualization of the distribution of the variable’s
values.
The density function (in blue in the example) is an
estimation of the probability distribution of the variable’s data. It’s
calculated by fitting different statistical distribution models to the
observed data. These models may include normal, Poisson, exponential
distributions, etc. The density function represents the probability that
an observation falls within a particular range of values.
The
probability distribution (in green in the example) represents the
statistical distribution model that would best match the density
function. It’s chosen by the user to best overlay with the blue curve.
The parameters of each probability distribution are approximated using
the mean and standard deviation of the variable.
By examining these three elements together, you can visually assess
how well the adjusted probability distribution model fits the observed
data. A good match between the three indicates that the model is a
precise representation of the data distribution. However, significant
discrepancies may indicate inadequacies in the chosen model or
particular characteristics of the data that require further analysis.
You can change the type of probability distribution to test which one
seems to fit best. If the chosen distribution does not fit at all, a
warning message appears. In the example, abundance is represented, and
the chosen law is a Lognormal distribution.
Once you are satisfied with the probability distribution, check the sentence above the “Modeling” button. There are two possibilities. If you have fewer than 30 observations, the sentence says: “You don’t have enough values to build a GLM or a GLMM” In this case, you need to change the working variable because there are not enough values to create a relevant model. Conversely, if the volume of data present in the dataset is sufficient (more than 30 observations), you will see: “Once you have chosen a distribution you can move on to building the model” When you are done, press the “Modelling” button.
Modelling
Model creation
This section is devoted to creating a model for inferential analysis.
The variable being analysed is shown at the top left of the tab and can
be modified in the “data representation” tab. The application allows you
to perform 3 types of inferential tests: GLMM, GLM and PERMANOVA. This
chapter begins with a reminder of the general statistical principles
used in the application.
General principles of statistics needed to use the tool
Parametric and non-parametric methods
Le domaine des statistiques existe parce qu’il est impossible de collecter des données auprès de tous les individus concernés (population). La seule solution consiste à collecter des données auprès d’un sous-ensemble (échantillon) des individus concernés, mais le véritable objectif est de connaître la “vérité” sur la population. La population est approchée en étudiant des variables descriptives. Chaque variable est un objet statistique, qui peut être décrit par des indicateurs. Les indicateurs statistiques tels que la moyenne, l’écart-type et les quartiles servent à résumer l’information concernant une variable observée. Lorsque l’on étudie un échantillon considéré comme représentatif, ces indicateurs servent à construire la loi de distribution de la variable étudiée. Chaque indicateur correspond à un “paramètre” de cette loi de distribution. On considère alors que la loi de distribution obtenue pour cette variable à partir de l’échantillon est applicable à la population. Étant donné que l’on ne peut généralement pas obtenir de données sur l’ensemble de la population, on ne peut pas connaître les valeurs des paramètres pour cette population. Il est toutefois possible de calculer des estimations de ces quantités pour l’échantillon. Lorsqu’elles sont calculées à partir des données de l’échantillon, ces quantités sont appelées “statistiques”. Une statistique estime un paramètre. Les procédures statistiques paramétriques reposent sur des hypothèses concernant la forme de la distribution (c’est-à-dire une distribution normale) dans la population sous-jacente et sur la forme ou les paramètres (c’est-à-dire les moyennes et les écarts-types) de la distribution supposée. Les procédures statistiques non paramétriques ne reposent sur aucune ou peu d’hypothèses concernant la forme ou les paramètres de la distribution de la variable dont l’échantillon a été tiré.
Linear models
Un modèle linéaire classique est une méthode paramétrique qui permet d’étudier la liaison statistique entre une variable réponse Y et les variables explicatives X. Soit yi la réponse de l’individu i et xi les valeurs prises par les variables explicatives pour cet individu. La relation entre X et Y peut s’écrire sous la forme : \[Y = α + \sum_jβ_jX_j + ε \] où ε représente les résidus du modèle, la variance de la variable Y non expliquée par les variables explicatives X, distribuée selon une loi normale d’espérance nulle. Le terme α correspond à ce qu’on appelle l’intercept et βj représente les coefficients estimés du modèle des variables explicatives Xj. La variable réponse pour un modèle linéaire doit être une variable approximativement normalement distribuée.
Generalized linear models
Les modèles linéaires trouvent une large application, mais ne peuvent
pas gérer des réponses continues clairement discrètes ou asymétriques.
Par exemple, les variables réponses de type “comptage”, souvent
asymétrique ainsi que les variables binaires comme la présence/absence
ne suivent pas une loi normale. Les modèles linéaires ne sont donc pas
adaptés à ce type de variables. Les modèles linéaires généralisés (GLM)
permettent l’extension des idées de modélisation linéaire à une classe
plus large de types de réponse, comme celles énoncées précédemment, sous
une méthodologie de modélisation commune. Une chose importante à
comprendre dans les GLM est la relation entre les valeurs de la variable
de réponse, Y (telles que mesurées dans les données et
prédites par le modèle dans les valeurs ajustées) et le prédicteur
linéaire. Le prédicteur linéaire émerge du modèle linéaire comme une
somme de chaque terme du modèle. Le prédicteur linéaire correspond à la
variable Y seulement lors d’un modèle linaire classique
suivant une loi normale. Dans le cas de modèle linéaire généralisé,
c’est la fonction de lien, g, qui relie la valeur
Y à son prédicteur linéaire N. \[ N =g(Y) \] La valeur de
N est obtenue en transformant la valeur de
Y par la fonction de liaison g, et la
valeur prédite de Y est obtenue en appliquant la
fonction de liaison inverse à N.
En utilisant
différentes lois de distribution et donc différentes fonctions de lien,
il est possible d’observer les conséquences sur les hypothèses des
résidus du modèle. La fonction de lien la plus
appropriée est celle qui produit les résidus les plus
conformes.
Generalized linear mixed model
Les modèles linéaires mixtes généralisés (GLMM) sont une extension des GLM. Un GLMM est dit “mixte” parce qu’il comprend au moins un effet “fixe”, les variables explicatives et au moins un effet “aléatoire”. Les effets aléatoires ne sont pas des termes évalués, ils servent uniquement à indiquer au modèle que les données ne sont pas indépendantes et reflètent une corrélation entre les unités statistiques. D’un point de vue statistique, cela permet d’estimer précisément la déviance résiduelle et donc d’éviter de biaiser l’erreur standard des paramètres. Au final, cela se traduit par des p-values plus fiables.
Analysis of variance using permutations
La PERMANOVA, ou Analyse de Variance Permutationale Multivariée, est une méthode statistique qui permet d’analyser les différences entre plusieurs groupes définis par des caractéristiques qualitatives, comme par exemple les différents traitements dans une étude expérimentale.
Contrairement à d’autres méthodes statistiques qui nécessitent certaines hypothèses sur la distribution des données, la PERMANOVA ne se base pas sur ces suppositions. Elle se focalise plutôt sur une matrice de distance entre les éléments étudiés. Cette approche lui permet de travailler avec des données de différentes dimensions, qu’elles soient simples ou complexes, et peu importe le nombre de catégories.
L’objectif de la PERMANOVA est de déterminer s’il existe des différences significatives dans la variabilité entre les groupes. Pour ce faire, elle évalue la variation entre les groupes (SS inter) par rapport à la variation à l’intérieur des groupes (SS intra). Une SS inter élevée suggère des différences importantes entre les moyennes des groupes, tandis qu’une SS intra faible indique une similarité accrue des observations au sein de chaque groupe.
La décision de rejeter ou non l’hypothèse nulle (l’absence de différence entre les groupes) se fait en comparant le rapport entre la variation inter-groupe et la variation intra-groupe à une distribution obtenue par permutation des données. Si ce rapport est significativement élevé, cela indique que les différences observées entre les groupes sont probablement réelles.
Cependant, la PERMANOVA présente quelques limites. Elle ne permet pas de déterminer quel groupe spécifique diffère des autres, seulement qu’au moins un groupe est différent. De plus, la présence de valeurs nulles peut biaiser l’estimation de la similarité entre les éléments, ce qui est particulièrement problématique en écologie où un zéro peut signifier l’absence d’une espèce. Cette limitation peut être atténuée en choisissant un coefficient d’association approprié dans le calcul de la matrice de distance.
Writing the model with GranulatShiny
Selon la méthode de modélisation que vous choisissez, la formulation du modèle diffère. Si l’on prend l’exemple de l’abondance comme variable réponse, le GLMM prendra en compte deux variables explicatives fixes, le traitement et la saison, et leur interaction ainsi que deux variables explicatives aléatoires, la campagne et la station : \[GLMM → Abun \sim traitement * saison + (1\|campagne) + (1\|station)\] Le GLM prendra en compte uniquement les variables explicatives fixes prises en compte dans le GLMM et leur interaction : \[GLM → Abun \sim traitement * saison\] La PERMANOVA prendra en compte les mêmes variables explicatives que le GLM : \[PERMANOVA → Abun \sim traitement * saison\]
Les modèles sont centrés sur la variable traitement car le suivi des
concessions d’extraction de granulats marins est basé sur la méthode
BACI (Before After Control Impact). Par définition, la méthode BACI
compare des sites témoins (c’est-à-dire non impactés) et des sites
impactés et teste les différences entre l’avant et l’après. Il s’agit
d’une méthode couramment utilisée dans la surveillance de
l’environnement océanique et une méthode BACI bien conçue reste l’une
des meilleures approches pour les programmes de surveillance des effets
sur l’environnement. Malheureusement, cette méthode présente plusieurs
limites qui compromettent sa capacité à détecter des effets notamment
parce que l’océan est spatialement et temporellement dynamique, et que
trouver deux emplacements statistiquement identiques l’un à l’autre tout
en étant suffisamment éloignés géographiquement pour être
statistiquement indépendants constitue un véritable défi.
Pour un GLMM et un GLM, vous devrez choisir une distribution de
probabilité. Par défaut, il propose la dernière distribution de
probabilité que vous avez vérifiée dans la partie précédente. Attention
la méthode utilisée pour la modélisation est itérative, il se peut donc
que la distribution qui semblait la plus adéquate dans la partie
précédente n’est pas forcément celle qui permettra de mieux faire
converger le modèle. Néanmoins l’onglet “diagnostic d’analyse” devrait
avoir permis de sélectionner un nombre de distribution possible pour ne
pas avoir à toutes les tester ici.
Vous pouvez également conserver ou non l’interaction entre les
covariables traitement et saison. Attention si l’interaction n’apporte
rien au modèle celle-ci est retirée automatiquement. Vous pouvez
également ajouter d’autres covariables dans votre modèle. Elles seront
ajoutées sans interaction avec les autres. Lorsque vous êtes prêt, vous
pouvez cliquer sur “démarrer la modélisation”.
Generalized linear mixed model
C’est la méthode à prioriser. Dans le paragraphe 8.4.2 du “Protocole
halieutique”, il est dit que pour évaluer la variabilité temporelle et
spatiale des différents indicateurs des ressources halieutiques avant
extraction, il faut utiliser des modèles linéaires généralisés à effets
mixtes (GLMMs) avec les variables temporelle et spatiale définies comme
effets aléatoires croisés plus un effet saisonnier fixe. La première
sortie est une reproduction de la sortie du logiciel r pour la ligne de
commande correspondante. Vous pouvez choisir d’afficher le tableau
d’analyse de la déviance, qui synthétise les résultats de la
modélisation en ne fournissant qu’un résumer de l’évaluation de
l’importance des effets fixes dans le modèle et de comprendre leur
impact sur la variable réponse. Vous pouvez également afficher le résumé
exhaustif de ces résultats. Vous pouvez choisir d’afficher les résultats
du modèle avant optimisation via le choix initial ou alors le modèle
optimisé via le choix final en bas à gauche.
Reproduction de la sortie R du modèle GLMM sur l’abondance.
## Generalized linear mixed model fit by maximum likelihood (Laplace
## Approximation) [glmerMod]
## Family: gaussian ( identity )
## Formula: log(Abun) ~ traitement * saison + (1 | campagne) + (1 | station)
## Data: dataset
##
## AIC BIC logLik deviance df.resid
## 1676.1 1726.4 -827.0 1654.1 708
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -3.1925 -0.6111 0.0535 0.5899 3.5320
##
## Random effects:
## Groups Name Variance Std.Dev.
## campagne (Intercept) 0.010794 0.10390
## station (Intercept) 0.004959 0.07042
## Residual 0.573266 0.75714
## Number of obs: 719, groups: campagne, 36; station, 20
##
## Fixed effects:
## Estimate Std. Error t value Pr(>|z|)
## (Intercept) 10.68015 0.10461 102.092 < 2e-16 ***
## traitementSans impact 0.11829 0.12315 0.961 0.33679
## saisonSpring 0.01607 0.14258 0.113 0.91025
## saisonSummer -0.15548 0.14264 -1.090 0.27571
## saisonAutumn -0.04808 0.14261 -0.337 0.73603
## traitementSans impact:saisonSpring 0.46811 0.16796 2.787 0.00532 **
## traitementSans impact:saisonSummer 1.05791 0.16815 6.291 3.14e-10 ***
## traitementSans impact:saisonAutumn 0.44990 0.16803 2.678 0.00742 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Correlation of Fixed Effects:
## (Intr) trtmSi ssnSpr ssnSmm ssnAtm trtmntSnsimpct:ssnSp
## trtmntSnsim -0.719
## saisonSprng -0.681 0.491
## saisonSummr -0.681 0.491 0.500
## saisonAutmn -0.681 0.490 0.500 0.499
## trtmntSnsimpct:ssnSp 0.491 -0.682 -0.720 -0.360 -0.360
## trtmntSnsimpct:ssnSm 0.491 -0.682 -0.359 -0.720 -0.359 0.499
## trtmntSim:A 0.490 -0.681 -0.360 -0.359 -0.720 0.500
## trtmntSnsimpct:ssnSm
## trtmntSnsim
## saisonSprng
## saisonSummr
## saisonAutmn
## trtmntSnsimpct:ssnSp
## trtmntSnsimpct:ssnSm
## trtmntSim:A 0.498La partie haute de la fenêtre de résultat rappelle le modèle qui a été utilisé pour calculer les effets.
## [[1]]
## [1] "Generalized linear mixed model fit by maximum likelihood (Laplace Approximation)"
##
## [[2]]
## [1] "gaussian" "identity"
##
## [[3]]
## glmer(formula = log(Abun) ~ traitement * saison + (1 | campagne) +
## (1 | station), data = dataset, family = gaussian(link = identity))Il y a le type de modèle avec la méthode de calcul utilisée. Après, il y a la loi de probabilité et sa fonction de lien. Enfin, la commande complete est affichée ce qui permet de vérifier que la bonne commande a été effectuée.
Ensuite, il est possible de lire des scores liés à la vraisemblance du modèle par rapport aux données et aux paramètres sélectionnés pour la construction du modèle.
Etant donné un échantillon observé (x1,…,xn) et une loi de probabilité Pθ, la vraisemblance quantifie la probabilité que les observations proviennent effectivement d’un échantillon (théorique) de la loi Pθ. On appelle vraisemblance associé à la loi de probabilité Pθ, la fonction L tel que : \[\displaystyle L(x_1,\ldots,x_n,\theta) = \prod_{i=1}^n P_\theta(x_i)\; \]
## AIC BIC logLik deviance df.resid
## 1676.0889 1726.4454 -827.0445 1654.0889 708.0000Ces indicateurs sont souvent fournis pour évaluer la qualité de
l’ajustement du modèle et aider à la sélection du meilleur modèle parmi
plusieurs candidats.
AIC (Critère d’information d’Akaike)
: L’AIC est un critère de sélection de modèle qui prend en
compte à la fois la qualité de l’ajustement du modèle et sa complexité.
Il favorise les modèles qui s’ajustent bien aux données tout en étant
simples. Un modèle avec un AIC plus bas est considéré comme préférable.
Cependant, l’AIC ne fournit pas d’indication sur l’ajustement absolu du
modèle, mais seulement sur son ajustement relatif par rapport aux autres
modèles candidats.
BIC (Critère d’information bayésien)
: Le BIC est un autre critère de sélection de modèle qui, comme
l’AIC, prend en compte à la fois l’ajustement et la complexité du
modèle. Cependant, le BIC pénalise plus sévèrement la complexité du
modèle que l’AIC. Un modèle avec un BIC plus bas est considéré comme
préférable. Contrairement à l’AIC, le BIC favorise la parcimonie, ce qui
signifie qu’il préfère les modèles plus simples.
logLik
(Log-vraisemblance) : La log-vraisemblance est une mesure de
l’ajustement du modèle aux données. Elle représente la probabilité que
les données observées soient générées par le modèle ajusté. Plus la
log-vraisemblance est élevée, meilleure est l’ajustement du modèle aux
données.
Deviance : La deviance est une mesure de
l’ajustement du modèle par rapport à un modèle de référence, souvent un
modèle nul. Elle est calculée comme la différence entre la déviance du
modèle ajusté et celle du modèle de référence. Une deviance plus faible
indique un meilleur ajustement du modèle aux données.
df.resid (degrés de liberté résiduels) : Les degrés de
liberté résiduels représentent le nombre de données indépendantes
restantes une fois que le modèle a été ajusté. Ils sont utilisés pour
calculer les statistiques de test et les valeurs p associées.
Les “scaled residuals” sont les résidus du modèle. Des tests sont effectués dessus afin de vérifier la bonne convergence et le bon ajustement du modèle.
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -3.192524 -0.611129 0.053534 0.000016 0.589887 3.532013Le tableau des effets aléaoires est spécifique au GLMM. Il renseigne
les informations sur cette partie de la formule :
(1 | campagne) +
(1 | station)
## Groups Name Std.Dev.
## campagne (Intercept) 0.10390
## station (Intercept) 0.07042
## Residual 0.75714Enfin, il y a la partie sur les effets fixes. Cette partie permet de dresser un diagnostic sur les facteurs et la variable d’étude. Le tableau des effets fixes fournit des informations clés sur les effets estimés des variables prédictives, leur précision et leur importance, aidant ainsi à comprendre les relations entre les variables et à tirer des conclusions sur les données.
Estimation (Estimate): Cette colonne indique les
coefficients estimés (ou effets) de chaque variable prédictive du
modèle. L’effet estimé de l’“Intercept” représente la valeur moyenne
estimée de la variable réponse lorsque toutes les autres variables
prédictives sont nulles.
Erreur standard (Std.
Error): Cette colonne indique les erreurs standard associées à
chaque estimation de coefficient. Les erreurs standard mesurent la
variabilité de l’estimation. Des erreurs standard plus faibles indiquent
des estimations plus précises.
Valeur t (t value):
Cette colonne indique la statistique t permettant de tester l’hypothèse
nulle selon laquelle le coefficient est égal à zéro. Elle est calculée
en divisant l’estimation par son erreur standard. Des valeurs t absolues
plus élevées indiquent une preuve plus forte contre l’hypothèse
nulle.
Pr(>|z|): Cette colonne indique la valeur
p associée à la statistique t pour chaque coefficient. Elle indique la
probabilité d’observer les données si l’hypothèse nulle (aucun effet)
était vraie. Des valeurs p plus faibles suggèrent une preuve plus forte
contre l’hypothèse nulle et indiquent que le coefficient est
statistiquement significatif.
## Estimate Std. Error t value
## (Intercept) 10.68014605 0.1046130 102.0919158
## traitementSans impact 0.11828840 0.1231490 0.9605307
## saisonSpring 0.01607218 0.1425775 0.1127259
## saisonSummer -0.15547873 0.1426408 -1.0900021
## saisonAutumn -0.04807517 0.1426082 -0.3371136
## traitementSans impact:saisonSpring 0.46810522 0.1679592 2.7870176
## traitementSans impact:saisonSummer 1.05790539 0.1681487 6.2914879
## traitementSans impact:saisonAutumn 0.44990449 0.1680289 2.6775420
## Pr(>|z|)
## (Intercept) 0.000000e+00
## traitementSans impact 3.367882e-01
## saisonSpring 9.102479e-01
## saisonSummer 2.757122e-01
## saisonAutumn 7.360312e-01
## traitementSans impact:saisonSpring 5.319559e-03
## traitementSans impact:saisonSummer 3.144374e-10
## traitementSans impact:saisonAutumn 7.416455e-03La colonne “Estimate” permet de déterminer la valeur moyenne que
prend la variable étudiée (dans notre exemple l’abondance totale) par
modalité des variables explicatives. Dans notre exemple, l’abondance est
expliquée par la variable “traitement” avec 2 modalités (impact, Sans
impact) et la variable “saison” avec 4 modalités (winter, spring,
summer, autumn). La ligne (Intercept) correspond à une valeur de base.
Cette valeur de base est associée à une modalité de chacune de nos
variables. Ainsi en hiver et avec une impact le logarithme de
l’abondance (car loi Lognormale) vaut en moyenne 10.68.
Pour
connaitre la valeur du logarithme de l’abondance en hiver et sans
impact, il faut additionner la valeur Intercept à la valeur de
l’estimate de la ligne “traitementSans impact” soit 10.68 + 0.12 qui
vaut 10.8. Et pour vérifier si le changement est significatif, il suffit
de regarder la colonne “Pr(>|z|)” et vérifier si la valeur est
inférieure à 0.05. Maintenant, pour obtenir la valeur du logarithme de
l’abondance en été et avec impact, il faut additionner la valeur
Intercept à la valeur de l’estimate de la ligne “saisonSummer” soit
10.68 + (-0.16) soit 10.52. Enfin, en additionnant la valeur Intercept
avec celle des lignes “traitementSans impact”, “saisonSummer”,
“traitementSans impact:saisonSummer”, il est possible d’obtenir la
valeur du logarithme de l’abondance en été sans impact soit 10.68 + 0.12
+ (-0.16) + 1.06 qui donne 11.70. Dans l’exemple, c’est l’interaction
entre la saison et le traitement qui apporte des changements
significativement différents. Ainsi, si ces deux variables étaient
considérées séparément, leur effet sur l’abondance ne serait pas
visible.
Et la dernière section montre les corrélations entre les termes à effet fixe du modèle. Chaque ligne et chaque colonne représentent un terme à effet fixe, et les valeurs du tableau sont les corrélations entre ces termes. Ces corrélations sont calculées sur la base de la matrice de covariance des estimations des effets fixes. Elles indiquent comment les effets estimés des différents facteurs fixes du modèle sont liés les uns aux autres.
Analyse a posteriori du modèle GLMM sur l’abondance.
Après avoir obtenu les résultats de la modélisation, l’application propose une série de graphiques pour diagnostiquer la qualité de l’ajustement du modèle, en se concentrant notamment sur l’analyse des résidus. Ces graphiques sont générés à l’aide du package DHARMa (“Residual Diagnostics for Hierarchical (Multi-level/Mixed) Regression Models”) dans R.
Le graphique de gauche, appelé “QQ plot residual”, est une représentation des résidus attendus par rapport aux observations réelles. Dans ce graphique, chaque point représente un résidu calculé par le modèle pour une observation donnée. Idéalement, ces points devraient suivre de près une ligne rouge diagonale, ce qui signifierait que les résidus sont distribués de manière approximativement normale. Si les points s’éloignent de manière significative de cette ligne rouge, cela suggère une mauvaise adéquation du modèle aux données observées.
En plus de la visualisation des résidus, l’outil DHARMa propose trois
tests pour évaluer la qualité de l’ajustement du modèle :
Test de Kolmogorov-Smirnov :
Ce test d’hypothèse
est utilisé pour évaluer si l’échantillon de résidus suit une loi de
distribution connue, déterminée par sa fonction de répartition continue.
Une déviation significative des résidus par rapport à cette distribution
attendue peut indiquer une inadéquation du modèle aux données.
Un test de Dispersion :
Ce test compare
l’écart-type observé des résidus à celui qui serait attendu en se basant
sur la simulation des données. Si la différence est significative, cela
peut suggérer une sous- ou sur-dispersion des résidus par rapport aux
attentes du modèle.
Un test de Valeur Aberrante :
Ce test vise à vérifier si le nombre d’observations dont les
résidus se trouvent en dehors de l’enveloppe de simulation est conforme
aux attentes du modèle. Une déviation significative de ce nombre peut
indiquer la présence de valeurs aberrantes ou une mauvaise adéquation du
modèle.
Chaque test fournit une mesure de la déviation par rapport aux attentes du modèle avec une p-value associée. Une p-value faible (< 0,05) indique généralement une déviation significative par rapport aux attentes du modèle, tandis qu’une p-value élevée suggère que la déviation observée pourrait être due au hasard et n’est pas statistiquement significative. Si cette déviation est significative, elle est signalée en rouge, indiquant que le test correspondant n’est pas conforme aux attentes du modèle. Ces diagnostics aident à identifier les inadéquations entre le modèle et les données observées et à guider les ajustements nécessaires pour obtenir un modèle plus approprié.
Sur le graphique de droite, des tests sont réalisés sur l’uniformité
et l’homogénéité de la variance des groupes évalués dans le modèle. Le
test “within-group deviation from uniformity” est un boxplot qui
représente la distribution des déviations des résidus au sein de chaque
groupe défini par les modalités des facteurs qualitatifs de votre
modèle. Chaque groupe est représenté par une boîte, où la médiane est
indiquée par une ligne à l’intérieur de la boîte, le premier et le
troisième quartile sont représentés par les bords inférieur et supérieur
de la boîte, et les moustaches s’étendent jusqu’aux valeurs maximale et
minimale. Les points au-delà de cette limite sont considérés comme des
valeurs aberrantes. L’objectif de ce test est d’identifier les groupes
pour lesquels les résidus présentent des variations importantes par
rapport à une distribution uniforme (ceux-ci apparaissent alors en
rouge). Des variations importantes peuvent indiquer une inadéquation du
modèle pour certains groupes spécifiques.
Le deuxième test
correspond à un test de Levene. Le test de Levene est utilisé pour
évaluer si les variances des résidus diffèrent significativement entre
les groupes définis par les modalités des facteurs qualitatifs. Il teste
l’hypothèse nulle selon laquelle les variances sont égales entre tous
les groupes. Une p-valeur faible (généralement < 0,05) indique une
différence significative dans les variances des résidus entre les
groupes, suggérant que l’hypothèse d’homogénéité des variances n’est pas
valide.
Dans le cas particulier où il y a trop de modalités différentes dû à de multiple covariables (ou que les variables explicatives soient quantitatives), la partie sur l’uniformité et l’homogénéité des groupes est remplacé par une représentation des résidus du modèle en comparaison avec les prédictions du modèle. Si n’y a aucun problème alors la phrase : “No significant problems detected” s’affiche en haut du graphique. Si des déviations des résidus par rapport à une distribution uniforme à travers différentes quantiles sont significatives ou que les déviations quantiles observées sont statistiquement significatives, les tests asociés apparaîssent en rouge, ce qui suggère une inadéquation du modèle pour certains aspects des données. Enfin, les valeurs aberrantes de la simulation (points de données qui se situent en dehors de la plage des valeurs simulées) sont mises en évidence par des étoiles rouges. Ces points doivent être interprétés avec précaution, car nous ne savons pas “à quel point” ces valeurs s’écartent des attentes du modèle. L’important est de vérifier que les tests de vérification ne soient pas significatifs.
Dans le cas de l’exemple, les tests de Kolmogorov-Smirnov, de
Valeur aberrante et de dispersion ne sont pas significatifs donc il n’y
a pas de problème. Si un de ces tests s’affichaient en rouge, cela
indiquerait que le modèle n’est pas optimal et il serait possible alors
de chercher un autre modèle qui s’ajusterait mieux. Comme ces modèles se
basent sur des données réelles, il est parfois impossible de trouver un
modèle parfait. Il faut alors choisir le modèle avec le moins
d’avertissement. On peut voir également que le test d’uniformité est
validé mais pas celui d’homogénéité. Une fois le modèle validé, vous
pouvez changer d’onglet et passer à la visualisation des effets associés
au modèle.
Si l’on rajoute la covariable année à notre modèle GLMM,
cela représente 96 modalités. Le graphique comparant les groupes est
remplacé par un graphique comparant les prédictions globales avec les
résidus du modèle.
Dans certains cas, la modélisation GLMM ne converge pas. Cela
signifie que les données disponibles ne permettent pas à l’algorithme de
calcul, associé à la formulation du modèle décidée par l’utilisateur,
d’estimer des valeurs de paramètres. Dans ces cas là, un message
d’erreur apparait : “Il y a une erreur lors de la modélisation.
Veuillez changer la loi ou le modèle.” Il est possible aussi que le
modèle produise des résultats mais dont l’analyse des résidus a
posteriori n’est pas satisfaisante.
Exemple :
Lorsque la modélisation par GLMM ne converge pas, il faut lui
préférer des méthodes de modélisation associées à des algorithmes de
calcul plus sipmles, c’est-à-dire avec moins de paramètres à estimer.
L’outil GranulatShiny en propose deux : le GLM et la PERMANOVA.
NOTA
BENE : lorsque le jeu de données utilisé pour la modélisation est
constitué de 30 observations non nulles ou moins, il est préférable de
s’en tenir à la méthode la moins coûteuse en termes de calcul, à savoir
la PERMANOVA.
Generalized linear model
Reproduction de la sortie R du modèle GLM sur l’abondance
##
## Call:
## glm(formula = log(Abun) ~ traitement * saison, family = gaussian(link = identity),
## data = dataset)
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 10.68271 0.09201 116.106 < 2e-16 ***
## traitementSans impact 0.11410 0.11770 0.969 0.33266
## saisonSpring 0.01560 0.13012 0.120 0.90458
## saisonSummer -0.14765 0.13012 -1.135 0.25687
## saisonAutumn -0.05388 0.13012 -0.414 0.67896
## traitementSans impact:saisonSpring 0.46888 0.16645 2.817 0.00498 **
## traitementSans impact:saisonSummer 1.04555 0.16660 6.276 6.04e-10 ***
## traitementSans impact:saisonAutumn 0.45938 0.16645 2.760 0.00593 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for gaussian family taken to be 0.5925863)
##
## Null deviance: 529.85 on 718 degrees of freedom
## Residual deviance: 421.33 on 711 degrees of freedom
## AIC: 1674.2
##
## Number of Fisher Scoring iterations: 2Le GLM sort des résultats proches de celui du GLMM et s’analyse de la même manière. Cependant, celui est moins précis. Il ne prend pas en compte les effets aléatoires induits par l’environnement ou la méthode utilisée.
Analysis of variance by permutation
Méthode de PERMANOVA intégrée dans l’application
Dans le cas de Granulatshiny, la PERMANOVA est appliquée sur une
matrice colonne regroupant des la plupart du temps des abondances ou des
biomasses. Ainsi les doubles zéros ne peuvent donc pas être pris en
compte dans le calcul de ressemblance. Donc les données sont
quantitatives et les doubles zéros ne sont pas considéré donc la méthode
qui semble la plus approprié est celle du coefficient de Bray-Curtis.
C’est le choix qui avait été fait à l’origine.
L’indice de
dissimilarité de Bray-Curtis, est utilisé en écologie et biologie pour
évaluer la dissimilarité entre deux échantillons donnés, en termes
d’abondance de taxons présents dans chacun de ces échantillons. Elle est
compris entre 0 (les deux échantillons ont la même composition) et 1
(les échantillons sont totalement dissemblables). La dissimilarité de
Bray-Curtis est souvent utilisée dans la littérature. Elle est
asymétrique et semimétrique. Elle se calcule comme ceci : \[ d_{jk}=\frac{\sum_{i} |x_{ij}-x_{ik}|}{\sum_{i}
(x_{ij}+x_{ik})} \] où i = colonne; j,k = lignes comparées; x
= valeurs d’abondances
Dans le cas où la matrice d’entrée n’a qu’une seule colonne soit une seule espèce et que dans les valeurs d’entrées, il existe des zéros, il arrive parfois que le dénominateur soit égale à zéro ce qui n’est pas possible et donc créer une erreur dans la matrice de distance. C’est le cas pour tous les indicateurs habituellement utilisés pour des données d’abondances qui pondèrent leur distance en fonction de l’abondance totale dans les sites comparés. Cette méthode ne pouvant s’utiliser dans notre cas, un autre coefficient de calcul de distance a été recherché. Celui-ci ne devait également pas tenir compte des doubles zéros sur des quantitatives.
La méthode retenue est la métrique du chi². Celle-ci donne davantage de poids aux espèces rares qu’aux espèces communes. Son utilisation est recommandée lorsque les espèces rares sont de bons indicateurs de conditions écologiques particulières. Pour appliquer cette méthode, il faut d’abord standardiser les données selon la méthode du chi² comme ceci : \[ x'_{ij}=\frac{x_{ij}}{\sum x_{j} * \sqrt \sum x_{i}} \] où i = colonne; j = ligne; x = valeurs d’abondances
Ensuite, on calcule la matrice de distance en calculant la distance euclidienne sur la matrice de données standardisées. \[d_{jk}=\sqrt \sum_{i} (x_{ij}-x_{ik})² \] où i = colonne; j,k = lignes comparées; x = valeurs d’abondances standardisées
L’inconvénient avec cette méthode, le calcul des effets se produisant sur des données transformées, il est n’est pas possible de quantifier directement l’impact d’un effet sur la variable initiale. On ne peut donc pas dire si un effet est plus ou moins fort sur la donnée initiale car celui-ci s’applique à la donnée transformée. Par contre, si un effet est considéré significatif sur les données transformées alors il l’est également sur les données initiales.
Reproduction de la sortie R de la PERMANOVA sur l’abondance
## Permutation test for adonis under reduced model
## Terms added sequentially (first to last)
## Permutation: free
## Number of permutations: 999
##
## adonis2(formula = dist ~ traitement * saison, data = dataset, permutations = 999)
## Df SumOfSqs R2 F Pr(>F)
## traitement 1 0.013749 0.08340 72.435 0.001 ***
## saison 3 0.008051 0.04883 14.138 0.001 ***
## traitement:saison 3 0.008098 0.04912 14.220 0.001 ***
## Residual 711 0.134958 0.81864
## Total 718 0.164856 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Dans le tableu de sortie de la PERMANOVA, on retrouve le nombre de
permutations et la formule du modèle. Ensuite on retouve plusieurs
indicateurs associés à chaque covariable explicative.
Df
(degrés de liberté): Cette colonne indique les degrés de
liberté associés à chaque terme du modèle.
SumOfSqs (Somme
des carrés): Cette colonne indique la somme des distances au
carré entre les observations dans l’espace multivarié.
R2
(R-carré): Cette colonne indique la proportion de variance
expliquée par chaque terme du modèle. Par exemple, pour “traitement”,
8,34 % de la variation des données peut être expliquée par le facteur
traitement.
F (statistique F) : Cette colonne
indique la statistique F pour chaque terme, qui vérifie si la variation
expliquée par ce terme est significativement plus importante que ce que
l’on attendrait du hasard. Des valeurs F plus élevées indiquent des
preuves plus solides contre l’hypothèse nulle d’absence d’effet.
Pr(>F) (valeur p): Cette colonne indique la valeur p
associée à la statistique F pour chaque terme. Elle indique la
probabilité d’observer les données si l’hypothèse nulle d’absence
d’effet (c’est-à-dire si toutes les moyennes des groupes sont égales)
était vraie. Des valeurs p plus faibles suggèrent une preuve plus forte
contre l’hypothèse nulle et indiquent que le terme est un prédicteur
significatif de la variation.
La PERMANOVA est la méthode à utiliser
en dernière. Elle esty moins précise et apporte moins d’informations
qu’un GLMM ou qu’un GLM.
Analyse complémentaire a posteriori de la PERMANOVA sur l’abondance.
Les mêmes boites à moustaches que dans l’onglet Représentation
des données sont affichées. Sauf que si la variable de comparaison (dans
l’exemple c’est le traitement) a un effet significatif sur la variable
expliquée (ici l’abondance) alors la p-value apparait en rouge en haut à
gauche du graphique. S’il n’y a pas d’effet détecté pendant la PERMANOVA
alors le message “Pas d’effet” apparait en haut à gauche du
graphique.
Effects representation
Cet onglet permet de visualiser graphiquement les effets des
variables explicatives sur la variable expliquée dans le cas d’un GLMM
ou d’un GLM. Dans le cas d’une PERMANOVA, cette section n’est pas
sollicitée et la fenêtre graphique sera blanche. Cette partie
retransforme les estimates du modèle en l’unité initiale (dans le cas de
l’abondance c’est un nombre par km²). Ainsi on peut voir la valeur
moyenne de l’abondance en fonction de la saison et du traitement.
D’abord vous devez choisir les deux prédicteurs à représenter.
Si vous avez plusieurs covariables vous devez les fixer afin de
pouvoir visualiser le graphique.
Dans l’exemple d’un GLM qui
regarde l’abundance totale en fonction du traitement et de la saison,
voici le graphique obtenu :
Ce graphique est une autre manière de représenter le tableau de
sortie du modèle de l’onglet précédent.
Statistical power
Cette partie est en cours de développement. L’outil antérieur
construit par Mathis Cambreling fonctionne seulement pour le jeu de
données ayant servi de base à ses calculs. L’outil n’étant pas
généralisable, celui-ci a été retiré pour assurer la stabilité actuelle
de l’application. Un autre outil est en cours de développement.
Bibliography
Anderson MJ (2017) Permutational Multivariate Analysis of Variance (PERMANOVA). Wiley StatsRef: Statistics Reference Online. John Wiley & Sons, Ltd, pp 1–15
Avezard C, Lavarde P, Pichon A, Legait B, Wallard I (2017) Impact environnemental et ́economique des activit ́es d’exploration ou d’exploitation des ressources minérales marines.
Bolker BM (2008) Ecological Models and Data in R. doi: 10.2307/j.ctvcm4g37
Bolker BM, Brooks ME, Clark CJ, Geange SW, Poulsen JR, Stevens MHH, White J-SS (2009) Generalized linear mixed models: a practical guide for ecology and evolution. Trends in Ecology & Evolution 24: 127–135
Colwell R (2009) Biodiversity: concepts, patterns, and measurement. The Princeton Guide to Ecology. pp 257–263
David V (2019) Statistique pour les sciences environnementales. ISTE Editions, Londres, Royaume-Uni
Gorodetska N, Behaghel G, Dalifard T, Daniel F, Grison X, Hausermann B, Laurent C, De Lantivy S, Lefebvre E, Panonacle H, et al (2023) L’ ́economie bleue en France.
Gregorius H-R, Gillet EM (2008) Generalized Simpson-diversity. Ecological Modelling 211: 90–96
Legendre P, Gallagher ED (2001) Ecologically meaningful transformations for ordination of species data. Oecologia 129: 271–280
Methratta ET (2020) Monitoring fisheries resources at offshore wind farms: BACI vs. BAG designs. ICES Journal of Marine Science 77: 890–900
Ministère de l’Environnement de l’́energie et de la mer (2016) Guide méthodologique pour l’élaboration des documents d’orientations pour une gestion durable des granulats marins (DOGGM). Ministère de l’Environnement, de l’Energie et de la Mer. Paris
MTE, UNPG, IFREMER, DREAL, DIRM (2023) Guide technique pour l’élaboration des ́etudes d’impact préalables à la recherche et l’exploitation des granulats marins. 48
Oksanen J (2022) Dissimilarity Indices for Community Ecologists.
Ortiz-Burgos S (2016) Shannon-Weaver Diversity Index. In MJ Kennish, ed, Encyclopedia of Estuaries. Springer Netherlands, Dordrecht, pp 572–573
Parent S-E (2020) Analyse et modélisation d’agroécosystèmes.
Rassweiler A, Okamoto DK, Reed DC, Kushner DJ, Schroeder DM, Lafferty KD (2021) Improving the ability of a BACI design to detect impacts within a kelp-forest community. Ecological Applications 31: e02304
Seger KD, Sousa-Lima R, Schmitter-Soto JJ, Urban ER (2021) Editorial: Before-After Control-Impact (BACI) Studies in the Ocean. Frontiers in Marine Science 8:
Shannon CE (1948) A mathematical theory of communication. The Bell System Technical Journal 27: 379–423
Smokorowski KE, Randall RG (2017) Cautions on using the Before-After-Control-Impact design in environmental effects monitoring programs. FACETS 2: 212–232
Underwood AJ (1994) On Beyond BACI: Sampling Designs that Might Reliably Detect Environmental Disturbances. Ecological Applications 4: 4–15
Walker R, Bokuniewicz H, Carlin D, Cato I, Dijkshoorn C, Backer AD, Dalfsen J van, Desprez M, Howe L, Robertsdottir BG, et al (2016) Effects of extraction of marine sediments on the marine environment 2005-2011. doi: 10.17895/ices.pub.5498
WGEXT (2019) Working Group on the Effects of Extraction of Marine Sediments on the Marine Ecosystem (WGEXT). doi: 10.17895/ices.pub.5